!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! calculate x4, y4 from given nozzle wall function
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! wall function is y=a+b*x+c*x^2
! slope of wall is dy/dx=tan(theta)=b+2*c*x
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! nomenclature:
!   a,b,c    wall equation constant
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! warning:
!   all angle in rad, not degree
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

subroutine Wall(lp)
   use VariableDef
   implicit none
   real*8::lp
!   integer::wallflag

! attachment point
   xa=rtd*sin(thetaa)
   ya=yt+rtd*(1-cos(thetaa))
   
   aw=((tan(thetaa)+tan(thetae))*xa**2-2*xa*xe*tan(thetaa)+2*ya*(xe-xa))/(2*(xe-xa))
   bw=(xe*tan(thetaa)-xa*tan(thetae))/(xe-xa)
   cw=(tan(thetae)-tan(thetaa))/(2*(xe-xa))
   
   ye=aw+bw*xe+cw*xe**2
   
   if (abs(thetae-thetaa) <= 0.01) then
      cw=0.0
      x4=(aw-y2+lp*x2)/(lp-bw)
      y4=aw+bw*x4+cw*x4**2
      theta4=atan(bw+2*cw*x4)
   else
      x4=(lp-bw-sqrt((lp-bw)**2-4*cw*(aw-y2+lp*x2)))/(2*cw)
      y4=aw+bw*x4+cw*x4**2
      theta4=atan(bw+2*cw*x4)
   end if
   
end subroutine Wall
